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Abstract 

We generalize the derivation of viscous anisotropic hydrodynamics from kinetic theory to allow for 
non-zero particle masses. The macroscopic theory is obtained by taking moments of the Boltzmann 
equation after expanding the distribution function around a spheroidally deformed local momentum 
distribution whose form has been generalized by the addition of a scalar field that accounts non- 
perturbatively (i.e. already at leading order) for bulk viscous effects. Hydrodynamic equations 
for the parameters of the leading-order distribution function and for the residual (next-to-leading 
order) dissipative flows are obtained from the three lowest moments of the Boltzmann equation. 
The approach is tested for a system undergoing (0+1 (-dimensional boost-invariant expansion for 
which the exact solution of the Boltzmann equation in relaxation time approximation is known. 
Nonconformal viscous anisotropic hydrodynamics is shown to approximate this exact solution more 
accurately than any other known hydrodynamic approximation. 

PACS numbers: 12.38.Mh, 25.75.-q, 24.10.Nz, 52.27.Ny, 51.10.+y 
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I. INTRODUCTION 


Relativistic fluid dynamics has been extensively used to describe the soft collective motion 
of relativistic heavy-ion collisions (see, for instance, and references therein) which plays 
a central role in the phenomenology of the quark-gluon plasma. This led to a number of 
works aimed at exploring the theoretical foundations of relativistic fluid dynamics, with the 
goal of identifying improved hydrodynamic approximations of the underlying microcopic dy¬ 
namics wm- For systems of deconfined quarks and gluons, the shear viscosity is expected 
to be much larger than the bulk viscosity at very high temperatures. It is then typically 
assumed that the bulk viscous pressure, whose Navier-Stokes value is porportional to the 
bulk viscosity, can be ignored in applications to heavy-ion physics. However, at temperature 
regimes typically produced experimentally in heavy-ion collisions, the order of magnitude 
of the bulk viscosity is unknown, and near the quark-hadron phase transition it could be 
large, due to the breaking of scale invariance by critical fluctuations and correlations [I3j [14] . 
Therefore, it is not well justified to a priori neglect the bulk viscous pressure when mod¬ 
eling the dynamics of QCD matter created at the Relativistic Heavy Ion Collider (RHIC) 
at Brookhaven National Laboratory and the Large Hadron Collider (LHC) at CERN. In¬ 
deed, several recent studies discussed the possibility of non-negligible bulk viscous effects 
on heavy-ion observables [T5H25] . One reason for possibly larger than originally expected 
bulk viscous effects in heavy-ion collisions is the existence of shear-bulk couplings in the 
equations of motion that control the evolution of the shear and bulk vicous pressures 
Heavy-ion collisions are characterized by initially very large differences between the longitu¬ 
dinal and transverse expansion rates that cause large shear stress which, in turn, generates 
bulk viscous pressure via bulk-shear coupling [26]. This mechanism should be included in 
phenomenological applications. 

A key assumption made in hydrodynamics is that the system remains close to local 
equilibrium. This assumption breaks down during the very early expansion stage of the 
systems formed in ultrarelativistic heavy-ion collisions. To account for these large early-time 
deviations from local momentum isotropy, a framework called “anisotropic hydrodynamics” 
(aHydro) was developed [23 EH] and recently generalized to “viscous (or second-order) 
anisotropic hydrodynamics” (vaHydro) [12] . In anisotropic hydrodynamics, one expands 
around an anisotropic background f a where the largest local momentum-space anisotropies 
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are built already into the leading-order (LO) term: 

f(x,p) = ,/a^ — + 5 f(x,p)- (1) 

Here S MI/ is a second-rank tensor whose structure depends on the shape and amount of the 
momentum-space anisotropy in the local fluid rest frame, and A is a temperature-like scale 
which reduces to the true local temperature in the isotropic equilibrium limit. Ref. 03 
was the first to include the correction 5f (which was computed using a Grad-Isreal-Stewart 
14-moment approximation), leading to the equations of viscous anisotropic hydrodynamics 
(vaHydro). When taking moments of the Boltzmann equation with the ansatz (Jl|, the 
contributions from 5f lead to additional dissipative (irreversible) currents n and that 
account for local momentum anisotropies not already built into the leading-order distribu¬ 
tion function f a (y/p-’E-p/A) and their relaxation equations |12j. In [12] these equations were 
eventually simplified by assuming massless degrees of freedom which allowed the longitudi¬ 
nal and transverse components Vl and V± of the anisotropic local pressure to factor into an 
isotropic thermal equilibrium pressure V lso multiplied by longitudinal and transverse “defor¬ 
mation factors” TZl,±{£) which depend only on the local momentum anisotropy parameter £. 
Of course, this assumption also implied zero bulk viscous pressure n. As a test of the vaHy¬ 
dro approach these simplified equations were then solved numerically for a transversally 
homogeneous system undergoing boost-invariant longitudinal expansion ((0+l)-dimensional 
expansion) for which the underlying Boltzmann equation can be solved exactly in the re¬ 
laxation time approximation (RTA). Comparison of this exact solution with vaHydro as 
well as several other viscous hydrodynamic approximations revealed a uniformly superior 
performance of the VAHYDRO scheme. 

We here generalize the vaHydro approach to the massive case where the above sim¬ 
plifications no longer hold. A suitable generalization of the tensor S^(x) in Eq. ([Tj) was 
first written down in [29] but not immediately exploited. It was recently shown that, with 
this generalization, already LO anisotropic hydrodynamics (without the 5f terms) implicitly 
contains some of the shear-bulk couplings present in modern versions of second-order viscous 
hydrodynamics [261 30]. By generalizing the work ra to the massive case, we here extend 
the works [26, [30] to next-to-leading order (he. we generalize the non-conformal aHydro 
approximation used in {26, 30] to non-conformal vaHydro), keeping all additional shear 
and bulk viscous corrections arising from the Sf term in Eq. ([!]). This improved approach is 
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again tested in a system undergoing (0+l)-dimensional boost-invariant expansion for which 
the exact solution of the RTA Boltzmann equation was recently extended to a gas of massive 
particles [ 51] . 

The structure of this paper is as follows. In Sec. [XT] we briefly review how to derive the 
macroscopic hydrodynamic variables from kinetic theory by expanding the local distribution 
function around a spheroidal momentum distribution. In Sec. Hi we write down the form 
of the leading order anisotropic background as well as the perturbations around this form. 


Then in Sec. IV we present the equations of motion for non-conformal viscous anisotropic 
hydrodynamics. In Sec. [V] we summarize the basic steps applied to obtain a closed form 
expression of the dynamical equations of motion for the residual dissipative currents, using 


the macroscopic conservation laws. Sec. VI we simplify these equations for the limiting 
case of a (0+l)-dimensional longitudinally boost-invariant expansion with finite particle 
masses and compare their numerical solution to the exact result from solving the microscopic 


Boltzmann equation. In Sec. VII we present our conclusions. 

Our notation is as follows: We use natural units h = ks = c = 1. The Minkowski metric 
tensor is g ,iu = diag(+, —, —, —). Greek indices run from 0 to 3 and Latin indices from 
1 to 3. The summation convention for repeated indices (Greek or Latin) is always used. 
Our tensor basis, in the local rest frame, is Xtf = = (1,0, 0,0), X[ l = x M = (0,1,0,0), 

X% = y 11 — (0, 0,1, 0), and Xg = z' J ‘ — (0, 0, 0,1). The transverse projection operator = 
-at at = <r —vAuA is used to project four-vectors and/or tensors into the space orthogonal 
to u O The notations A ^ | and A^ = | (A tJjU —A utl ) denote symmetrization 

and antisymmetrization, respectively. A^ = A^A a ^ where = A^A^ — A^A Qi g/3 
is the transverse (to u ) and traceless projector for second-rank tensors. The four-derivative 
is = d/dx 11 , D = is the convective derivative (the time derivative in the comoving 
frame), V A ‘ = A^ u d u is the covariant notation for the spatial gradient operator in the local 
rest frame, and 6 = d is the scalar expansion rate. 


II. HYDRODYNAMICS FROM RELATIVISTIC KINETIC THEORY 

To keep the presentation selfcontained, we here briefly review how to extract hydrody¬ 
namic variables from the Boltzmann equation for an expansion around a locally spheroidal 
momentum distribution m In kinetic theory the one-particle distribution function / is 
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governed by the Boltzmann equation, 


P“dj = C[f] , ( 2 ) 

where C[f] is the collision kernel. The particle current and energy-momentum tensor are 
expressed as the first and second moments of the one-particle distribution function 

J' 1 = <?/*) , , (3) 

where we defined the average of a momentum-dependent observable 0(p) at point x as 

(0)(x) = J dP 0(p)f(x,p) (4) 

with the Lorentz invariant momentum-space measure dP = (2n)~ 3 (d 3 p/E). We decompose 
the particle four-momentum p M into parts parallel and orthogonal to the four-velocity vP of 
the local fluid rest frame [32], 

= Eu> x + Pi X», (5) 

where E = p-u(x) is the local rest frame energy and pi= — X^{x)p v are the Cartesian com¬ 
ponents of the four-momentum in the local rest frame. For systems that are locally approx¬ 
imately spheroidal in momentum-space, characterized by a “local anisotropic equilibrium” 
distribution function / a , we can decompose / as 

f(x,p) = U{x,p) +Sf . (6) 

Then and can be tensor decomposed as 

J ■** = A/V* + , (7) 

= £u»u v - {V± + fl)A^ + (V L - V±) z^z u + . (8) 

Here z^{x) = X£{x) is the four-vector that reduces in the local fluid rest frame to a unit 
vector in longitudinal (z) direction, J\f is the particle density and V^ is the particle current 
in the local rest frame, £ is the energy density in the local rest frame, V± and Vl are the 
transverse and longitudinal pressures, n is the bulk viscous pressure, and is the shear 
stress tensor defined by 

U = <£)., V“ = (p.) s AT, 

£={E 2 )„ P± = (p1>. = (9) 

n = _i(A«Vp s ) s , ^ = (pV>), . 
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In these equations we introduced the notation (• • • ) a = JdP(---)f a and {■■■)$ = 
f dP (■ ■ ■ )Sf, and made use of the generalized Landau matching conditions (E)$ = ( E 2 )^ = 
0. For later convenience, the total bulk viscous pressure II is calculable as 

n = 2P± 3 +Pl - v eq + n, (io) 


and the total shear stress tensor is obtained from 


Av v 

^ = (V l -V_l) ( — + z»z v ) + = (V L -V L ) 


x^x u +y^y u —2z >l z v 


+ 7f^. (11) 


III. 14-MOMENT APPROXIMATION FOR THE DISTRIBUTION FUNCTION 


A. Leading order (LO) distribution function 


In this paper we consider systems that are, to leading order, spheroidal in momentum- 
space in the local rest frame. This is accomplished by introducing the anisotropy tensor 
so that the leading-order one-particle distribution function takes the form: 


/» 



( 12 ) 


where we have assumed zero chemical potential and fo has the functional form of a local 
thermal equilibrium distribution, 


fo(y) = 


(13) 


ey + a ’ 

where a = ±1,0 corresponds to Fermi-Dirac, Bose-Einstein, and classical Boltzmann statis¬ 
tics, respectively. 

The most general decomposition of the rank-two tensor zL ,lu that possesses spheroidal 
symmetry in the local rest frame is 


5^ = vPu v - ± iz^z v , 


(14) 


where all terms are functions of position x. In local rest frame coordinates / a takes the form 

/a = fo Q \J m 2 ±(l±$)pi±(l±4>±0P^ = fo , (15) 

where we defined E 2 = (l+$)mj_ cosh 2 y ± ^m 2 ± sinli 2 y — $m 2 , with m 2 L = rnf+pf. Con¬ 
necting the “anisotropic equilibrium” quantities with moments of / a , one hnds that these 
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quantities can be written as (see Appendix) 


AT (A ,f,$;m) = <£> a = 

£ (A,£,<&;m) = (£ 2 ) a 
A 4 


Ao(A) 


( i +$)( i +$+0 ’ 


( 16 ) 


= - f dy cosh 2 y f dm± m\ fo(^yJrti±[( l+$) cosh 2 y + £, sink 2 y] — , (17) 


7MA,£,*;m) = <pl> 
A 4 


= —j j dy j drh± rh± m 2 ) [( 1 +$) cosh 2 y + £sinh 2 y] — ?h 2< hj, (18) 


(A,f,$;m) = (pjj) a 

^4 /-oo 
* '0 

where rh = m /A and m_i_ = m±/A= ^Jm 2 -\-p\/ h. The equilibrium thermodynamic quanti¬ 
ties are given as moments of /q. For Boltzmann statistics they have the functional form: 


= -^- 5 - f ciy sinh 2 y f drh± m\ f 0 (^yjm 2 ± [(!+<!>) cosh 2 y + £ sinh 2 j/j — ?h 2 $ j, (19) 


' T ' 777 

Ao(T; m) = —K 2 (m/T ) , 


T 2 m 2 / 771 \ 

£ 0 (T;m) = — r (3AT 2 (m/T) + -Ad(m/T)) , 


2 vr 2 

P 0 (T;m) =AA(T;m)T, 

where K n (z) are the modified Bessel function of the second kind. 

B. 14-moment expansion of the deviation from the LO distribution 


( 20 ) 

( 21 ) 

( 22 ) 


In the 14-moment approximation, the deviation 5f of the full distribution function / 


from the locally anisotropic state ( 12 ) is expanded to second order in momenta as m 


= a - /3 E + wE 2 - , 

Ja/a ^ 


(23) 


where / a = l— af a . In the absence of a chemical potential, as assumed in (12), there is no 
heat current V and the coefficients of any terms linear in in Eq. (23) vanish. By 


inserting Eq. (23) into the definitions ([9]) of the residual dissipative flows, the 14-moment 
coefficients can be expressed in terms of these flows by solving the matrix equation Ab = c, 
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where 


r Ji,o 

O 

1 

Jsfl + j'i,\ 

0 

0 

0 

J2,0 

—Jsfi 

J\$ + Ja,\ 

0 

0 

0 

J2A 

~Js, 1 

Ja,\ + §<^4,2 

0 

0 

0 

P X 2 X 1 

-<pl! 

Ptl + ¥>I 

A 1111 

0 

0 

0 

0 

0 

0 

2A 1212 

0 

0 

0 

0 

0 

0 

2A 1313 

¥>21 

-P 31 

Pll + ¥>42 

A 1122 

0 

0 

0 

0 

0 

0 

0 

0 

\P Z 2 1 

-pH 

¥>41 + ¥>42 

A 1133 

0 

0 





\ T 


a j3 

W Wn 

W \2 W 13 VJ-22 

W 23 

w 33 ) 



(O 0 II 7fll 7Tl2 7Tl3 ^22 ^23 ^33 


0 0 

0 0 

0 0 

A 1122 0 

0 0 

0 0 

A 1111 0 

0 2A 1313 

A 1133 0 


Pi) 

Plfl 

P2,l 

A 1133 

0 

0 

A 1133 

0 

\ 3333 


(24) 


(25) 

(26) 


This allows the distribution function expanded around an anisotropic background to be 
expressed in terms of the residual dissipative flows II and 7 R" as 


/ — /a + 


A n n + A£%„+ (Arn + Ar a ^, 


'n 


a/3 ) P(^Pu) 




(27) 


where An, AA^, and Aalong with the auxiliary tensors and A^""^ appearing 


in (24), are defined in Ref. 


IV. VISCOUS ANISOTROPIC HYDRODYNAMIC EQUATIONS OF MOTION 


In this section we derive the hydrodynamic equations of motion by taking moments of 
the Boltzmann equation. Taking moments implies multiplying (j2| by integer powers of the 
four-momentum and integrating over momentum-space. This process results in the following 
n-tli (n > 0) moment equation: 

d^ip 111 ■ ■ ■ p^ +1 ) = . (28) 


The n-th rank collisional tensor is defined in the following manner: 

= [ dPE r pP 1 ■ --p^Clf ], 


(29) 


with The inhnite set of coupled moments (28) is equivalent to knowing the 

full solution / of the Boltzmann equation (J2| . Only the first few moments have an intuitive 








physical meaning: The zeroth moment d fi (p fl ) = C embodies the conservation of particle 
number for vanishing C, the first moment d fl (p fl p u ) = 0 the conservation of energy and 
momentum. The macroscopic equations of spheroidal viscous anisotropic hydrodynamics 
are derived in the following subsections. 

A. Zeroth moment of the Boltzmann equation 

The zeroth moment of the Boltzmann equation gives 

d tl f = DU + U9 + dU' 1 = C. 

Denoting the action of the time derivative in the local rest frame D by a dot, Eq 
be written as an equation of motion for the rest frame particle density U: 

u = -ue - d,y ,j ' + c. 


(30) 


(30) can 


(31) 


B. First moment of the Boltzmann equation 


The first moment of the Boltzmann equation is equivalent to the requirement of energy- 
momentum conservation: d^T^ 11 = 0. With the viscous anisotropic hydrodynamic decom¬ 
position of given in ([7]) this conservation law yields 

d,jT llu = u’'D(E+V ± +n) + u v (E+V±+fl)9 + (E+V±+fl)Du v - d v (V±+fl) 

(32) 

+ z u D l {V l -V ± ) + z u (V l ~V ± )9 l + ( V L -V ± )D L z I/ + d^ v = 0. 

Projecting these four equations on the fluid four-velocity yields an equation of motion for 
the rest frame energy density £\ 

u v d t U ,v = £ + (£+V ± +fl)9 + (V L -V ± )u u D z z" - U v a lw = 0. (33) 


The projections A a ll d ll T ,lu transverse to u ^ yield equations of motion for the fluid four- 
velocity tdh 


= ( s+v ± +fi)u a - v a 0Pu+n) + ( 34 ) 

+z a D z (V L -V ± ) + z a (VL~V ± )(d^) + (V L -V ± )D z z a - {V^-Vi_)u a u v D z z v = 0. 
In the above equations we have introduced the velocity shear tensor a^ and the 


derivative operator Di = X^d^. Equations (33) and (34) are the fundamental equations of 
relativistic viscous anisotropic hydrodynamics. 
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C. Second moment of the Boltzmann equation 


The second moment of the Boltzmann equation gives 


d,,r wX = C vX , 


( 35 ) 


where J r ^ uX = (p^p v p x ). Decomposing p M into parts parallel and orthogonal to u M by using 
Eq. ([5]) leads to 

T» vX = (E?)vPu v u x + (E 2 Pi ) {u^u u X x + v?X\u x + X?u v u x ) 

+ (Epipj) (u*XYX x + X'lu v X x + X?Xy) + (p iPjPk )X?X"X x . (36) 


To evaluate the l.li.s. of Eq. (35) requires taking the four-divergence of tensor J z ^ vX 


d fl E lwX = D(E 3 )u u u x + (E 3 ) ( u v u x Q + 2 u (v Du x) ) 

+ X\X: x D(E PlPj ) + (EpiPj) (X''X X 6 + 2 x[ v DXf) 

+ 2 u^X^DiiEpiPj) + 2 (E PiPj ) + u^DiX?) . 


(37) 

Projecting out the transverse to tA and traceless part of Eq. (35), yields 

(D(£p,p,) + (EpiPj)&) + 2{Ep,p,)Af, (X'M* + X'DX) = C<“» . (38) 


To work out the averages ( Ep^j ) appearing in Eq. (38), we use (27) to write 


(E PtPj ) = ijj + vsn + , 


(39) 


where Z(q is defined in the Appendix and 




dP XnEpiPjfJ a + / dPp^p^EpiPjfJz 


n, 


dP Ep ^/ a / a + A 


Oipfj,!; 


dPp {a pp)Ep iPj fJ a 


7r 




(40) 

(41) 


Using (39), equation (38) can finally be written in the following form: 


xfxf 


IS + A 3 n + + sgn + 


+ vfxf 


ii; + t/’gn + *»■'#, 


/LiU 1 


iS + ffin + J a “f (.Y.-B.y + = c<“«. 


(42) 
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V. EQUATIONS OF MOTION FOR THE RESIDUAL DISSIPATIVE FLOWS 


To close the system of equations we need to know the space-time evolution of the dissi¬ 
pative currents appearing in ([£]). To accomplish this, we derive the equations of motion for 
II and 7f M " from their kinetic definitions mum-- 

,2 r 


m 


n = —— J dPSf, 

= f dPp (y0) S f . 


(43) 

(44) 


The Boltzmann equation (J2]) can be written in the form 


5/ = -/a--(p-V(/ a +5/)-C[/] 


(45) 


Substituting this into the expressions (43) and (44) one obtains the following equations of 
motion: 


—ri - c_i = w - xe - + —fie - 


m 


m* 


fill V 


(46) 


~(^) _ c (y) = + c »v + + n ^x ^ + u «v xZa ) + (! + YQ^ Xa ^xu a 

- -7f ^9 - 2vfJV " >a + 2 jr^u v)x + 

3 

- (E- 2 p^p u) P {a) P if)) )s V a up ■ (47) 


The evolution equations for the dissipative flows II and can now be obtained by inserting 


the closed form (27) of the single-particle distribution function into the expectation values 


(• • •) j on the r.h.s. of the equations of motion (46) and (47). This leads to 


n-c_i = w 

m z 


xe - y^a, u + X-m - Cunv,/^ _ (48) 


_ c (y) = ^u + c „u + + W vX ^ + /Va4 ) + (]_ + <1 >)Q^ Xa V X U a 

- -Tt^e - 2 Tt^(T v)x + 2 TT^u v)x + 21 la^ 

3 

- n- KT^xV a up . (49) 
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The dissipative forces for the bulk viscous pressure and shear-stress tensor evolution equation 
are defined as 


W = /3 a i7o,0,l + i 


X = — 

3 L 


y ,w = 


M iMl = ^<f> tel,-! AT x; + |A 


z^z v 


(50) 

(51) 

(52) 

(53) 


The remaining dissipative forces in Eq. (49) /C M ", C^, etc. and transport coefficients 


etc. are given in Appendix C of Ref. 

The viscous anisotropic hydrodynamic framework for a nonconformal system (with van¬ 


ishing chemical potential) in a general (3+l)-dimensional framework is defined by Eqs. (31), 
(33), (34), (42), (|48|), and (|49|). Structurally they reduce to Eqs. (88) and (90) of Ref. [12] 


when taking the limit $ —> 0. The difference between the equations studied here and in 
is that here we account for some of the bulk viscous effects non-perturbatively, by including 
them via the scalar field $ already in the LO distribution function / a . This leads to slight 
changes in the structure of the relaxation equations for II and and also changes the 
values of the transport coefficients. 


VI. (0+1)-DIMENSIONAL EXPANSION FOR A NONCONFORMAL SYSTEM 
A. Reduced evolution equations 

In this section we present and solve the boost-invariant VAHYDRO equations for a sim¬ 
plified situation without transverse expansion. In the following we will use the relaxation 
time approximation (RTA) for the scattering kernel, 


C{f] = [/(p;A,«,*)-/ 0 («-p/T)], 

7~eq 


(54) 


where r eq is the relaxation time, assumed to be momentum-independent. For transversely 
homogeneous systems undergoing boost-invariant longitudinal expansion, the Boltzmann 
equation (J2]) with an RTA collision kernel (54) can be solved exactly [334435] . and this can 
be used to determine the efficacy of various approximation schemes. In the situation just 
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described there are no transverse derivatives, the comoving time derivative A = DA simply 
becomes dA/dr, and the shear stress tensor is fully defined by a single non-vanishing 
component fr = 7 t z z — — 7 i zz \ at z — 0, n ,l v = diag(0, — 7 t/2 , — ff/2,7r). For (0+l)-dimensional 
expansion with azimuthal symmetry we have the following simplifications: 


x\D x u x = v\D y u x = 0 , 
z\D-u x - -t , 0 - - . 


(55) 

(56) 


Using this we can write the (0+l)-d viscous anisotropic hydrodynamic equations of motion 
as: 


Af=--~— (N-Neq) , 


T T, 


eq 


(57) 

(58) 


£ = — (y + Vl + n — y 

^{EvD - A&A) = ; (( E A) - 3(Ep 2 J) + A ((Ep 2 J - (Ep 2 )) . (59) 


dr 


T 


'eq 


With some algebra, using the explicit functional form (20) of the particle density (16) for 
Boltzmann statistics, the first two of these can be rewritten in terms of the parameters of 
the leading-order distribution f a as 

£ 


!+$+£ 


_ 2 I 3 ! m ^1 ( m /A) \ A / 2 | 1 A * 

+ A K 2 (m/A) J A + Vl+$ + l +$+{J 


{d^£)i + (d A £)A + (<%,£)$ = -~(£ + V L + U-n) , 


• (60) 
(61) 


where the partial derivatives of £ on the left hand side of Eq. (61) can be worked out from 
the explicit expression The effective temperature T in Eq. ( J6Q| ) is obtained from the 

dynamical Landau matching condition 


£(A,f,4>; m) = £ e ,(T; m). 


(62) 


Some additional work yields the evolution equations for If and 7 r in the form 


n = -r 


(2V ± +V L _ , ^ , m 2 ( ~ A 1 ~ 22 • 3 ~ A , l+d>+£ 




-Peq+hl ] +— ^ Jo,0,1 ^ - 2^0A-J _ g^ 2 ’ 1 - 1< ^ > + 


rrzz 


-A --A * 

/ 'nn / 'ii7r ? 

r r 


(63) 
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7T = -r ( 7T - ~(V L ~ V±) 


1 

A 


Ji 


zz 

0 , 0,1 


*4+ 


rrzzzz _ rrzz 

.^0,0,-l ^2,1,-1 


3 ( j-zz 5 ~ 

- I 1+2,1,—! q ‘-'4,2,—1 


II 7T 

"t~ ^Trll T ^7T7T * 

r r 


$ 


(64) 


Equations (59)-(64) form the coupled set of dynamical equations that must be solved for 


(0+l)-dimensional expansion. The J integrals appearing in the last two of these equations 


are defined in the Appendix. The terms (Ep 2 ) and (. Ep 2 Z ) appearing in Eq. (59) involve 
the transport coefficients -0^ and ip% zz . From a formal point of view, it is nice to have 
analytic expressions for the transport coefficients. However, for nonconformal systems the 
“shear-bulk” coupling is rather complicated. For numerical purposes it is then easier to just 


use the parametrization of the non-equilibrium distribution function (23) and numerically 


invert the matrix equation Ab = c at each time step in the numerical integration to obtain 
the coefficients in the 14-moment approximation. We now show how to do this. For (0+1)- 
dimensional expansion with azimuthal symmetry we have the following simplifications: 


b = (a P w wu 0 0 u >22 0 w 33 ) , 


= 0 0 n jf/2 0 0 tt/2 0 - 7T 


x\D x u x = y\D y u x = 0 , 
Z\D z u x = -- , 9 — - 

T T 


By defining the vector 


/ T n = ( 'Tii 'Tii 'Tii _j_ T’ii 'jiixx p\ p\ T’iixx p\ T’iizz \ 

u ~ W 1,0 J2,0 ^ 3,0 ' '-'3,1 ^1,0 u u ^1,0 u ^1,0 J ’ 


we can write 


(Ep‘) = I“„ + b - J" . 


(65) 

( 66 ) 

(67) 

( 68 ) 

(69) 

(70) 


Then taking the derivative 


d o, dl {)o ( ,_ x dA,\ _,,dn 


Tr^ = 2TM _ j + + W+, (71) 


where we introduced the vectors 


u= 001000000 


V=(000|00i0-1 


(72) 

(73) 
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With the anisotropic form (12) as the underlying LO distribution function, it is convenient 


to evolve the system in terms of the kinematical parameters £, $, and A, rather than the 


macroscopic densities. Writing Eq. (71) in terms of £, and A 


h (Eft) = + v-Sft + . 


(74) 


where we have introduced the shorthand notation: 


« = a a % - b ) j“ + b ■ d a j “, 

V>n = (-4- 1 *) • J “. 

>/'“ = (.4 '„) . J7- fi . 


(75) 

(76) 

(77) 


Then dehning ipk = V’aT — V’aT; where fc G {£, A, <f>, II, 7r}, and taking the ^-component of 


Eq. (59), we get 


V’ff + V'aA + + V’nn + ip= - (( Ep 2 x ) - 3(Epl)) + — ({Epl) - ( Ep\ )) . (78) 


r 


r, 




The transport coefficients appearing in Eqs. (48) and (49) can be written as 

,2 


m 

T 


4 

Avi-n — g + 1 

where &n = .A -1 !!, b n = A~ l v , and 


o 

m 

(79) 

bnAJ n , 'W = —b^AT n , 

7 

' 3 7r 5 ^'K'K g + 7r ? 

(80) 


77n 

^TT 


'7’ZZ /TZZ I ^TZZ 'TXXZZ fl fl ‘TXXZZ fl ‘TZZZZ \ 

<-4-2,0; <-4-1,0; >-40,0 “r <-4o,l; <-4-2,0 ; u ; u ; >-4-2,0 ; u ; <-4-2,0 J ’ 

„ ~ e; 

/TZZ_ rrzzzz rrZZ _ rrzzzz /7ZZj__ rrzz _ rrzzzz _ rrzzzz 

<-4o,l <-4-2,0 ; <-4i,i <-4-1,0 ; <-42, 1 -1- g <-42,2 <-4o,0 <^ 0,1 ; 

rrXXZZ _ rrXXZZZZ n n rrXXZZ _ rrXXZZZZ ri rrzzzz _ rrZZZZZZ\ 

<-4o, 1 — *-4-2,0 ; U, U, *-4 q,i J_ 2,0 ; u ; «-4o,l <-4-2,0 J 


(81) 


(82) 


B. Numerical results 


In this subsection we solve for (0+l)-dimensional expansion the vaHydro equations 


(59)-(64) numerically and compare the resulting evolution histories for the macroscopic 


thermodynamic quantities with the corresponding moments of the related exact solution of 
the Boltzmann equation [31] . We initialize the system at To = 0.5fm/c with To = 600 MeV, 
II 0 = 0, and 7 f 0 = O.For simplicity and illustration, we assume a temperature independent 
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f 0 =0, 7'o=0.6 GeV, m=l GeV, r eq =0.5 fm/c f 0 =0, T 0 =0.6 GeV, m = l GeV, r eq =0.5 fm/c 




fo=100, T 0 = 0.6 GeV, m=l GeV, r eq =0.5 fm/c 


fo = 100, T o =0.6 GeV, m=l GeV, r eq =0.5 fm/c 




FIG. 1: Ratio of the longitudinal to transverse pressure V^/Vl (left column) and the bulk viscous 
pressure II (right column). The top panels correspond to an initial anisotropy parameter £o = 0 
whereas the bottom panels are for £o = 100. The black solid, red short-dashed, blue dashed-dotted, 
and green long-dashed lines are the results obtained from the exact solution of the Boltzmann 
equation, NLO anisotropic hydrodynamics (vaHydro), LO anisotropic hydrodynamics (aHydro), 
and second-order viscous hydrodynamics [Ml, respectively. The initial conditions in this figure 
are To = 600 MeV, m = 1 GeV, flo = 0, ttq = 0, T eq = 0.5 fm/c, and To = 0.5 fm/c. 

relaxation time, exploring the cases T eq = 0.5 fm/c in Figs. [I] and [2] as well as the ten times 
larger value r eq = 5 fm/c in Fig. [3] 

In Figs. M we plot in the left panels the evolution of the ratio V±/Vl between the 
longitudinal and transverse pressures, and in the right panels the bulk viscous pressure II. 
For Figs, [l] and [3] we assume particles of mass m — 1 GeV while in Fig. [2]we use m — 0.1 GeV. 
The upper and lower panels in Figs, [l] and [2] correspond to different initial momentum-space 
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f 0 =0, T 0 = 0.6 GeV, m =0.1 GeV, r eq =0.5 frn/c f o =0, T 0 = 0.6 GeV, m=0.1 GeV, r eq =0.5 fm/c 




fo=100, 7 0 = 0.6 GeV, m=0.1 GeV, r eq =0.5 fm/c f 0 = 100, T 0 =0.6 GeV, m =0.1 GeV, r eq =0.5 fm/c 




FIG. 2: Similar to Fig. [2j but for a ten times smaller mass m = 0.1 GeV. 

anisotropies: For the histories shown in the upper panels of Figs. [Tj and [2j as well as those 
in Fig. [3j we assumed initial momentum isotropy, £o = 0. The lower panels in Figs. [l]and 
[2] start instead from a very anisotropic initial state with £ 0 — 100. In all cases, the solid 
black line shows the results obtained from the exact solution of the Boltzmann equation. The 
short-dashed red lines represent our VAHYDRO results while the dashed-dotted blue and long- 
dashed green curves show results from LO anisotropic hydrodynamics (aHydro) [30] and 
second-order viscous hydrodynamics (vHydro) evaluated in the 14-moment approximation 
M for comparison. 

We observe that in all cases the VAHYDRO framework gives the best approximation 
to the exact solution for the longitudinal/transverse pressure ratio and the bulk viscous 
pressure. Especially during the early evolution stages, vaHydro matches the exact solution 
almost perfectly while LO aHydro and second-order viscous hydrodynamics [7H9] exhibit 
significant deviations. The improvement achieved by including in the dynamical description 
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£o=0, T 0 =0.6 GeV, m=l GeV, r eq = 5 fm/c £ 0 = 0, r 0 =0.6 GeV, m = 1 GeV, r eq =5 fm/c 




FIG. 3: Similar to Fig. [IJ but for r eq = 5 fm/c. An additional set of purple dotted curves shows the 
effect of setting the bulk-shear coupling terms in the evolution equations for II and 7 r to zero. 


the additional dissipative flows generated by Sf is particularly visible in the evolution of the 
bulk viscous pressure, caused by the non-vanishing particle mass. 

As discussed in [26J, for massive particle systems coupling terms between the shear and 
bulk viscous pressures play an important role in the evolution of the viscous stress; in viscous 
hydrodynamics (vHydro) these must be included explicitly at second order in an expansion 
around a locally isotropic momentum distribution |7H9] while anisotropic hydrodynamics 


(aHydro), based on the ansatz (12), captures their effects already at leading order, with 
similar precision. The dotted purple curves in Fig. [3j which were obtained by setting in 


Eqs. (63) and (64) the bulk-shear coupling coefficients An^ and A^n to zero, show that in 
our improved VAHYDRO framework residual bulk-shear coupling terms between II and 7r 
(due to 5f) play only a minor role, and only at late times. The main improvement over 
AHYDRO and second-order VHYDRO results from the other terms on the right hand sides of 
these equations, including the diagonal couplings Ann and \ nn . 


VII. CONCLUSIONS 

In this paper we derived a generalization of the viscous anisotropic hydrodynamic frame¬ 
work pjj] to systems with massive degrees of freedom, assuming a vanishing chemical poten¬ 
tial. To test the efficacy of this extended vaHydro formalism we applied it to a transversally 
homogeneous non-conformal system that undergoes boost-invariant longitudinal expansion 
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(i.e. (0+l)-dimensional flow) for which there exists an exact solution of the RTA Boltzmann 
equation mums]. We tested its precision in this situation by comparing, over a wide range 
of particle masses, relaxation times and initial momentum anisotropy parameters, the nu¬ 
merical predictions of VAHYDRO for the longitudinal/transverse pressure ratio and the bulk 
viscous pressure with the corresponding results obtained from the exact (0+l)-d solution of 
the RTA Boltzmann equation EHEB], as well as with those from two other hydrodynamic 
expansion schemes, namely second-order viscous hydrodynamics in the 14-moment approx¬ 
imation (vHydro) [7H3] and anisotropic hydrodynamics (aHydro) [HU] . In all cases we 
found that VAHYDRO agrees almost perfectly with the exact kinetic solution and presents a 
significant improvement over the other two hydrodynamic approaches. 

For massless theories, a powerful test for the efficiency of various hydrodynamic ap¬ 
proximation schemes for the RTA Boltzmann equation is the amount of entropy generated 
by viscous heating during the evolution of the system [I2l [29j. Among all known hydrody¬ 
namic schemes, only aHydro and VAHYDRO reproduce this quantity qualitatively correctly 
in both the strong (r re i —>• 0) and weak (r re i —>• oo) coupling limits, and for (0+l)-d expan¬ 
sion VAHYDRO does so almost perfectly, ft would have been nice to perform this check also 
for the case of massive theories. However, for massive particles the entropy needs to be 
computed from kinetic theory (rather than using the particle density). We postpone this for 
future work. 
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Appendix: Evaluation of the “anisotropic” thermodynamic integrals 


In this Appendix we compute the auxiliary “anisotropic” thermodynamic integrals 
and Jn qr ' lk defined as 

y - * = t2 ^+ 1) l)!! / dP El Pil ---p ilPh ---p jr U, (83) 

= (2 ( ~ + 1) 1 ' |, f APE ’’- 2, (A«^) 5 ■ -p jr fj. , (84) 

where i = 1, 2 and j = 3 denote the number of p* and p 2 components, respectively, and 
/ a = 1—a/a- To work out the derivative of I nqr Hn Jp with respect to £, $ and A we use the 
following relations: 


^ (£ a 7a) 

3* (Kh) 

d\f a 


2 E a \E a Aj /a/a ’ 


A al3 p a pfj f _ 1\ 

2 E a \E a Ay 





This allows us to write 


= Z I r 2_ rrh--Hjl-jp+2 

L/ <i- L 'nqr r> \ 1 *^n,q,r— 2 a ‘-'r 


A 


' n,q,r— 1 


\ • • • 1 -• • • • 
_ _ I -T _ /rti—iejl—jp 

{J Q- L 'nqr 2 l ^ n+2,q,r—2 ^n+2,<?,r*-l 


(85) 

( 86 ) 

(87) 

( 88 ) 

(89) 


fii-W-ip _ _J_ Pfii -Hh 'ip 
u A-~ L 'nqr j^2 ^n,q,r+l 

We parametrize the four-momenta in hyperbolic coordinates, 

y/ 1 = (mi cosh p, p_|_ cos (j), p_i sin 0, mi sinh p) , 


(90) 


(91) 


where m 2 = m 2 + p(^_ is the transverse mass, with integration measure 

dy m±dm±d(j) 
d = (2^)3 ' 


Then 


(92) 


= A-T— —t— j dym±dm±d<l>(m±cmhy) n ^ 2q (m 2 —mi cosh 2 j/)’ 
x hi ((l+$) cosh 2 p + £ sinh 2 y) — m 2< f>] (rn\— m 2 )^ 2 ■ ■ ■ p i( (m± sinh y) p 

x /o f ^ \/ m i ((1+3*) cosh 2 p + £ sinh 2 p) — m 2< lA , (93) 
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where we defined the scaled transverse Cartesian momentum components pi = Pi/p±, such 
that p x = cos (p and p y = sin0. We now define dimensionless parameters rh± = rn±/ A and 
m = m /A which results in 

( — 1)1 \n+t+p+r+2 r 

Inqr H3l "' 3P = (2q + 1 )!! - 2V 1 2 - ^ ^ J ^ ^ V V ' Hn qtpr ^ ^ ’ ( 94 ) 

where 

^ V " ie = J o ( 95 ) 

roc 

'H nq e pr (u, £, rh) = / dm 1 _m n ^ 2q+p+l \fh\ ((l+<3>) cosh 2 y + £ sinh 2 y) — m 2 4 $] r 2 (96) 

J rh 

x (m 2 — m 2 cosh 2 7 8 y) 9 10 11 (m 2 — m 2 )^ 2 f 0 (^\J™ ± ((1+$) cosh 2 y + £ sinh 2 y) — . 

The same decomposition follows for Jn qr H3l '" 3v \ 

( — 1)1 \n+£+p+r+2 r 

Jnqr" H3l " 3P = ( 2q + 1)\\ -^^ ^ J ^ COSh " _2<? V sinhP V ^nqipr (p, f m) , (97) 

where 'H nq ( pr takes into account quantum statistics and is obtained by making the substi¬ 
tution /o(-) —> /o(-)/o( - )- We note that in the classical limit (a = 0) the two functions are 
identical, 0£ itil "' jp (a=0 ) = J^' idl '" ip (a=0 ). 
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